function x=constrain(x,xl,xu)
% constrains a number x in [- infty , + infty] to [xl,xu]

x=(xu-xl)/(4*atan(1))*atan(x)+(xl+xu)/2;